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Abstract:  This  report  presents  a  summary  introduction  to  the  homogeni¬ 
zation  procedure  in  numerical  methods  via  sequential  projections  onto 
nested  subspaces  spanned  by  mutually  orthogonal  scaling  and  wavelet 
orthonormal  families  of  functions.  The  ideas  behind  the  technique  of 
multi-resolution  analysis  unfold  from  the  theory  of  linear  operators  in 
Hilbert  spaces.  The  homogenization  procedure  through  successive  multi¬ 
resolution  projections  is  presented,  followed  by  a  numerical  example  of 
sequential  analysis  and  synthesis  of  a  simple  signal  illustrating  the  appli¬ 
cation  of  the  theory.  A  structural  example  shows  a  practical  application  of 
multi-resolution  analysis  to  the  displacement  response  of  a  cantilever  with 
highly  heterogeneous  elasticity  subjected  to  a  concentrated  load  at  the  tip. 
An  introductory  appendix  describes  the  reproducing  kernel  methods  of 
mathematical  representation  of  a  given  field. 


DISCLAIMER:  The  contents  of  this  report  are  not  to  be  used  for  advertising,  publication,  or  promotional  purposes. 
Citation  of  trade  names  does  not  constitute  an  official  endorsement  or  approval  of  the  use  of  such  commercial  products. 
All  product  names  and  trademarks  cited  are  the  property  of  their  respective  owners.  The  findings  of  this  report  are  not  to 
be  construed  as  an  official  Department  of  the  Army  position  unless  so  designated  by  other  authorized  documents. 
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1  Introduction 

The  objective  of  this  report  is  to  provide  a  summary  introduction  to  the 
homogenization  procedure  in  numerical  methods  via  sequential  projec¬ 
tions  onto  nested  subspaces  spanned  by  mutually  orthogonal  scaling  and 
wavelet  orthonormal  families  of  functions.  The  development  is  intended  to 
be  essentially  self-contained.  The  mathematical  (Greenberg  1978;  Gilbert 
2006)  and  signal  processing  (Strang  and  Nguyen  1995;  Vetterli  and 
Kovacevic  1995)  backgrounds  are  outlined  in  logical  sequence  avoiding 
statement  and  theorem  demonstrations.  The  organization  of  the  mathe¬ 
matical  background  is  built  after  the  work  of  Gilbert.  Readers  needing 
further  expansions  of  the  exposed  material  should  consult  the  funda¬ 
mental  literature  listed  in  the  references. 

The  ideas  behind  multi-resolution  analysis  unfold  from  the  theory  of  linear 
operators  in  Hilbert  spaces  (Davis  1975;  Kolmogorov  and  Fomin  1975; 
Akhiezer  and  Glazman  1993;  Zwillinger  1996).  The  homogenization  pro¬ 
cedure  through  successive  multi-resolution  projections  is  presented  in 
Chapter  5  (Mehraeen  and  Chen  2004),  followed  by  a  numerical  example  of 
sequential  analysis  and  synthesis  of  a  simple  signal  illustrating  the  appli¬ 
cation  of  the  theory  (Wavelet  extension  of  Mathcad  2000).  Then,  a  struc¬ 
tural  example  shows  a  practical  application  of  multi-resolution  analysis  to 
the  displacement  response  of  a  cantilever  with  highly  heterogeneous 
elasticity  subjected  to  a  concentrated  load  at  the  tip.  The  report  ends  with 
a  succinct  introductory  appendix  describing  the  reproducing  kernel 
methods  of  mathematical  representation  of  a  subject  field  under  investi¬ 
gation  (Liu  and  Chen  1995). 
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2  Mathematical  Background 

Energy 

Let  Z  denote  the  set  of  integers:  Z=  {...,-2,-i,o,+i,+2,...}, 

and  R  denote  the  collection  of  real  numbers:  R  = 

The  energy  of  an  analogue  signal//)  defined  on  X  =  R  is  an  average 
measure  of  its  size  given  by  the  square  of  its  norm: 

E(f)=\\f(x)f=  p| /Ml  dx  (la) 

ii  ii  u —oo  '  ' 

in  a  one-dimensional  space,  or  by 

E(f)=  JJD\f{x,y)\2  dxdy  (lb) 


when/(v)  is  defined  in  a  two-dimensional  space  D  =  R2  (case  of  a  static 
image). 

Similarly,  the  energy  of  a  digital  signal  a  (i.e.,  a  sequence)  defined  onl  =  Z 
(or  on  D  =  Z2  in  the  two-dimensional  case)  is  given  by 

+00 

E(a)=  W  (2a) 

n=— oo 

or  by 

+00 

£(a)=  ^  amn  |  (2b) 

m,n=— oo 

in  one-  and  two-dimensional  spaces,  respectively. 

Then,  for  a  finite-energy  function  (or  signal)  one  has  and,  similarly, 

for  a  finite-energy  sequence  one  has  E(a)<°°.  In  the  analogue  case,  the  set 
of  all  finite-energy  functions  (i.e.,  functions  that  are  square-integrable)  is 
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denoted  by  LfX);  whereas,  in  the  digital  case,  the  set  of  all  finite-energy 
sequences  (i.e.,  sequences  that  are  square  summable)  is  denoted  by  L(Z). 

Inner  product  spaces 

The  energy  of  individual  elements  of  a  function  vector  space  can  be  con¬ 
veniently  calculated  using  the  notion  of  inner  product. 

Definition.  Let  v  be  a  function  vector  space.  An  inner  product  on  v  is  a 
mapping  that  translates  pairs  of  elements  of  the  space  into  corresponding 


scalars,  in  such  a  way  that 

(/,/)> 0,  and  (/,/)=  0  iff  f  =  0  (3a) 

{f,9)={9>f)>  (3b) 

(f  +  g,h)=(f,h)  +  (g,h)  (3c) 

(X-f,g)  =  X-(f,g )  =>  (f,X-g)  =  ( A •  g,f )  =  X-(f,g)  (3d) 


for  all/,  g,  and  h  in  v,  and  where  an  over-lined  variable  means  the  complex 
conjugate  of  the  variable.  Spaces  of  finite-energy  signals  become  inner- 
product  spaces  under  the  appropriate  inner-product  definition.  For  exam¬ 
ple,  for  the  pair  of  analogue  signals  f  g  in  L2{X)  one  has 

(/,£/)  =  f+  f{x)-g{x)dx  (4a) 

— oo 

and  for  the  pair  of  sequences,  a  =  (an}n  and  b  =  {bn}n  in  l2(Z),  one  has 

+OQ  _ 

(a,b)=  an\  (4b) 

71= — OO 

The  three  main  properties  of  a  complex  inner  product  in  v  are 
Cauchy-Schwarz  inequality: 

\{f,gf  <E(f)-E(g) 


(5a) 
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Triangle  inequality: 

Ey*(f+g)<E1/i(f)  +  EL/*(g)  (5b) 

Polarization  identity: 

(f,g)  =  —  { [E(f+g)  -  E{f-g)\  + *•  {E(f+i-g)  -E(f-i-g)] }  (5c) 

Orthonormal  families 

The  inner  product  if,  g )  of  two  functions  gives  a  measure  of  how  well  f 
correlates  with  g.  In  fact,  when/looks  like  g,  then  if,  g)  will  be  close  to  the 
energy  of  g. 

Definition.  A  subset  {cpn}  of  an  inner-product  space  v  is  an  orthonormal 
family  when 


(cp  i,<pj)=8ij  (6) 

where  8 y  is  Kronecker’s  delta,  or,  equivalently,  when  the  family  consists  of 
mutually  orthogonal  elements,  each  having  unit  energy.  For  example,  con¬ 
sider  the  family  of  sequences  {e^}n  ,  with  elements 

e(n)=(...  0,0, 1,0,0,...)  (7) 

with  unity  at  the  nth  location,  then 

(e«\eu>)  =  8s  (8) 


and  any  sequence  (a  =  (an}n)  e  ifZ)  may  be  represented  as 

a  =  £>n.s(n)  =  £(a>e(">y">  (9) 

n  n 

Definition.  The  expansion  of  a  function/ in  an  inner-product  space  v  in 
terms  of  the  orthonormal  family  {(pn}  in  v  is  defined  by 
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Slf\  =  Z)(/>(Pn)(Pn  (10) 

n 


The  numerical  values  (f,  cp„)  are  called  the  coefficients  of/. 

Notice  that  the  family  {cpn}  is  not  necessarily  complete,  i.e.,  a  basis  span¬ 
ning  the  space  v. 

Bessel’s  inequality 

Let  {cpn}  be  an  orthonormal  family  in  an  inner-product  space  v.  Then,  the 
inequality 


£|(/,<P„)f<-E(/)  (11) 

77 

holds  for  all/in  v.  In  other  words,  the  sequence  {  (f,  cpn)  }n  of  the  coeffi¬ 
cients  of /has  a  finite  energy  that  is  bounded  by  E(f)  and,  therefore, 
belongs  to  f2CZ). 

Definition.  An  orthonormal  family  {cpn}  in  v  is  complete  when 

£|(/.<P»)f  =  E(f)  (12) 


for  all/in  v. 

Given  a  complete  orthonormal  family  {cpn},  then  one  can  expand 

s{f]  =  J2(f’(Vn)<Vn=f  (13) 


in  the  mean-square  sense.  In  other  words,  the  distance  D(f,  S[/])  tends  to 
zero  in  the  limit: 


lim  E 


N 


f~  E  (/>(Pi)(P; 


i=—N 


—  0 ,  asN 


oo 


(14) 


Theorem.  If  (cpn}is  a  complete  orthonormal  family  in  an  inner-product 
space  v,  then  the  coefficients  (/cpn)  satisfy  the  identities: 
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( Plancherel )  E{f)  =  J2  {f>Vn 


(■ Parseval )  (f,g)  =  ^(/,  cpj-  (g,  cpn 


Fundamental  bases 

1.  The  simplest  basis  is  the  standard  basis  for  U{Zn)'. 

£(n)  =  (0,0,..., 0,1,0,..., 0,0) ,  —  oo  <  n  <  +oo 

a  complete  orthogonal  family  whose  nth  point  carries  the  only  non¬ 
component  (unity)  at  the  nth  location. 

2.  The  Fourier  basis  for  f 2(Zn): 


:(*) 


Vat 


a  k  2k  k-(N-l) 

1  ,  CO  ,  CO  ,...,  CO  1  j 


k  >0 


where 


•  231 
1\  N) 

co  =  e  {  ’  =  cis 


2n 


N 


is  an  N^1  root  of  unity,  in  the  sense  that  cofc  JV  =  l.  Thus, 


:(0) 


Vat 


:(D_  1  (1  rw,2  .JJV-1) 

D  ^  ?  0),  0)  ,...,0) 


:(N- 1)  _  1 

"  ~  yfN 


l^-Sco2^-15,...,^-^ 


(15a) 

(15b) 

(16) 

zero 

(17) 

(18) 

(19a) 

(19b) 

(19c) 


is  a  complete  orthonormal  family. 
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3.  Defining  the  fundamental  functions: 

(n)  _  1 

V2’ 

then  each  of  the  sets  { cp(n) }  and  { cp(n) }  is  orthonormal  in  (2  and,  together, 
the  set  { cp(n) ,  cp('!) }  is  a  complete  orthonormal  family  in  f2.  Consequently, 
every  point  x  in  f2  admits  the  block  expansion: 

x  =  0  W"»)  <p‘”’  +  ^(^cp00)  <pw  (21) 

77  n 

4.  Since 

fa  +  cis(2nm£)-cis(2nnl)d£  =  8mn  (22) 

the  family  {  cis(2jin^);  -oo<n<+oo  }  of  period  one  is  orthonormal  and  com¬ 
plete  for  the  space  L2[a,  a+i],  for  any  choice  of  a. 

5.  Since 


Jin 


,2n+l 


(20) 


cp(,!)  = 


e2  n  +  &2n+1 


and  cp 


/+°°  /  \  - 7 — C 

b0n(x)\n(x)dx  =  8mn  (23) 

-00 

the  family  {bon:-°°<n<+n},  bon  =  b(x-ri)  of  all  integer  delays  of  the  box 
function  with  scale  zero  is  orthonormal  in  L2[-°o,+oo],  However,  the  family 
is  not  complete.  (Consider  the  expansion  of  a  Haar  function  in  terms  of 
this  family,  for  example.) 

6.  Similarly,  the  Haar  family 

hkn{x)  =  2^  h(2k-x  —  nj,  -oo<n<+oo  (24) 

is  orthonormal,  but  not  complete,  in  L2[-°o,+oo]5  for  each  fixed  resolution  k. 
(Consider  the  expansion  of  a  box  function  in  terms  of  the  family,  for 
example.) 

7.  By  allowing  all  possible  resolutions  and  integer  delays  of  the  Haar 
function,  one  obtains  a  complete  orthonormal  family  in  L2[-°o,+oo]5  since 
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CM  „ 

2  2  ■/. 


+  00 


h\2J  ■x-m\- h(2  -x-l)  dx  =  8jk  dml 


(25) 


Every  finite-energy  function  in  L2[-°°,+°o]  admits  a  series  expansion: 

/(x)=S(/’MMx)  (26) 

k,l 


Conservative  mappings 

A  linear  operator,  3 :  v  =>  W,  from  one  inner-product  space  v  into  a  pos¬ 
sibly  different  inner-product  space  W  is  conservative  when  the  squared 
norm  of  every  point  in  the  space  is  preserved  upon  transformation. 
Conservative  transformations  also  preserve  inner  products. 
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3  Signal  Processing  Background 

Delay 

The  delay  operator  on  f2  is  defined  as  the  mapping 

^{^n-lK  (27) 

which  preserves  energy,  i.e., 

£(Wx)  =  2|x„_1|2=^|x„|2=i;(x)  (28) 

n  n 

Convolution 

The  convolution  operator,  or  filtering,  between  two  sequences  is  defined 
by  the  sequence  with  elements 

(  h  *X)n  =  XX  • X(n-m )  =  XX \n-m) '  Xm  (29) 

771  771 

or,  in  operator  terms,  as 

=  ?""(*)  (30) 

777 

By  the  triangle  inequality,  one  has 

E%{h*X)=E%  (  y,k  (3i> 

,771  )  [  777  J 

which  is  bounded  if  j  ^  hm  | j  =  K( a  real  constant).  Therefore,  the 

convolution  operation  preserves  energy  when  the  coefficients  of  the 
sequence  h  are  absolutely  convergent. 
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To  identify  the  adjoint  of  the  operator  h,  we  notice  that 


(h  *  X ,  y)  =  E  K-m) '  Xm  '  Vn  =  Yj  Xm  '  I]  \n-m)  '  Vn 


=  (x,h**y  (32) 


i.e.,  the  adjoint  operator  is  given  by  the  sequence  h*  =  j/i_n  j  . 

Discrete-Time  Fourier  Transform 

The  Discrete-Time  Fourier  Transform  (DTFT)  is  defined  by  the  pair  of 
inner  products: 


5e 


1  1 
2 ’2 


(33) 


xn=jVy2X(iyelln]dl 


(34) 


where  e\  =  {  cis(2nriQ  }n  .  Notice  that  the  transform  X(£),  also  called  the 
frequency  response  function  of  x,  has  period  one. 

The  DTFT  (A)  as  a  mapping  from  f2  into  L2[-i/2,  +1/2]  is  energy  preserv¬ 
ing,  and  maps  delay  into  modulation,  i.e., 

(v  xf  (I)  =  £>„_1  ■(T2"'”6  =  X(l)  ■  c!'s(2jt|)  (35) 


A  more  crucial  property  of  the  DTFT  is  that  it  maps  convolution  into 
pointwise  multiplication,  i.e., 


{h*x)A(l)  =  J2  \J2hn-m-Xm  '  ds(2jlTll)  =  H  $)  ■  X  (%)  (36) 


The  frequency  response  function  of  the  adjoint  h*  is  given  by 


h)  &)  =  \Y,h-n-cis(27inty\  =  Hfe) 


(37) 
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Since  the  transformation  is  conservative,  inner  products  remain  invariant 
upon  transformation,  i.e., 

(■ x>y )  =  /  y x(9  • Y ©  d ?  =  E[  f  yx ©  •  «s(2™0  •  yn 

A  _  "[  A  j  (38) 

=XX-i/n 

rz 

for  all  finite-energy  sequences  in  f2.  Consequently,  one  has 

xn  =  f^(0'«s(2Jin5)c^  (39) 

z-Transform 

The  z-transform  of  a  digital  signal  x  =  { xn}n  is  defined  as 

X(z)  =  ^xn-z“n  ,  zg  Complex  (40) 

n 

Since  on  the  unit  circle  in  the  complex  plane  z  =  z0  =  • cis  (2ji^),  one  can 
write 

xA  (?)  =  (*, e?)  =  •  «s(-2nn?)  =  J>n  •  z~n  (41) 

n  n 

and  the  z-transform  can  be  regarded  as  the  extension  of  the  DTFT(x)  from 
the  unit  circle  to  the  whole  complex  plane. 

Down-sampling 

The  down-sampling  operation  on  a  digital  signal  x  =  (x„}n  produces  a  new 
signal 

(|2)x  =  {x2„}n  (42) 

by  discarding  all  odd-indexed  terms  from  x  and  re-indexing.  Clearly,  the 
down-sampling  operation  in  not  energy  preserving,  since 

£ ( (l  2)  x)  =  £ |  x2n  |2  <  £ |  x„  |2  =  £  ( x )  (43) 

n  n 
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In  Fourier  transform  terms  one  has 

2'{  X  \  +X  ^1}  =  ~{Exn-cJs(n-^)  +  ^xn-ds(n-^  +  1))  + 

=  ~|  Exn-cl's(«-^)-[1  +  (-1)n]j  (44) 

and,  consequently,  one  obtains 


which  is  an  average  in  the  frequency  domain. 

Up-sampling 

By  contrast,  up-sampling  is  the  converse  of  down-sampling.  Given  a  digital 
signal  y  =  {yn}n,  up-sampling  produces  a  new  signal,  defined  by 

(T2 )y  =  {vn}n,  such  that  (v2n^yn)and(v2n+1^-  0)  (46) 

Indeed,  up-sampling  is  the  adjoint  operation  to  down-sampling,  since 

((1 2)  x,  y)  =  Ex2n  •  yn  =  (x,  (T  2)  y)  (47) 

77 

Up-sampling  is  a  conservative  operation  and,  in  terms  of  the  Fourier 
transform,  one  has 

(t  2)*  y]  (I)  =  ((t  2)y,  e5)  =  (y,  e(20)  =  Y( 21)  (48) 

which  is  a  V2-period  function  in  the  frequency  domain. 
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Filters 

A  discrete  filter  is  a  linear,  time-invariant  operator  E  acting  on  f2.  In  other 
words,  the  filtering  operator  E  satisfies  the  commutative  property,  i.e., 

E(¥i)  =  ¥(E(i))  (49) 

with  respect  to  the  delay  operator.  Any  such  operator  is  given  the  following 
convolution: 


h  *  x  =  J  (50) 

in  which  h  =  {hn}n  is  the  sequence  of  filter  coefficients.  Another  way  of 
expressing  this  crucial  operation  is 

<51> 

e 

When  only  finitely  many  of  the  ht  are  different  from  zero,  the  operator  E  is 
said  to  be  of  finite  impulse  response  (FIR). 

Since  the  DTFT  maps  convolution  into  point-wise  multiplication,  one  has 

y  —  E(x)  =  h*  x  =>  =  (52) 

where  H(£)  =  (h,  ei)  is  the  frequency  response  function  of  the  filter  coeffi¬ 
cients.  The  effect  of  E  is  to  select  or  to  reject  certain  frequencies  in  a  signal. 
An  FIR  filter  is  said  to  “try  hard”  to  be  a  low-pass  filter  when  H{ o)  ^  o  and 
H(±V 2)  =  o;  by  contrast,  E  is  said  to  “try  hard”  to  be  a  high-pass  filter 
when  H( o)  =  o  and  H(±V 2)  ^  o.  Corresponding  to  every  low-pass  filter  E, 
there  is  an  associate  high-pass  filter  E .  Define  E  as  the  FIR  filter  having 
coefficients  hn ,  such  that 


(53) 


Then,  by  Fourier  transforming  to  the  frequency  domain,  one  has 
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=  V  vy^1/ 


Some  fundamental  filter  examples  are 
l.  Haar’s  filters,  defined  by 


X/5’  "  =  ° 
K=  W^,  n  =  1 

0 ,  n  ^  0, 1 


(55a) 


^72  ’  71  =  0 

■/i'  "=1 

0  ,  n  ^  0, 1 


(55b) 


The  corresponding  frequency  response  functions  are 


H(l)  =  y[: 2  •  cz's(ji^)  •  cos(jri;) 


(56a) 


#(?)  =  V2 -z- czs (ji£) -sin (ji£). 


(56b) 


Clearly,  H  (Q  tries  hard  to  be  a  high -pass  filter  having  a  single  zero  at  the 
origin  in  the  frequency  domain. 

2.  The  tent  filters,  defined  by 


V  r,  n  =  0,2 
/2-V2 


<?n  =  B 


n  =  1 


A/2’ 

0 ,  otherwise 


(51a) 
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and 


9n  = 


,  n  =  — 1,1 

71  =  0 
otherwise 


The  corresponding  frequency  response  functions  are 

G©  =  V 2  •  cis (2j©  •  cos2  (j© 

and 

G©  =  V 2  •  sin2  (j© 


(51b) 


(58a) 


(58b) 


Clearly,  G©  tries  hard  to  be  a  high-pass  filter  having  a  double  zero  at  the 
origin  in  the  frequency  domain. 

(3)  Daubechies’s  db-2  filter,  defined  by  the  following  only  non-zero 
coefficients, 


1  +  V3 
!°  4V2  ’ 


3  +  V3 


4V2  ’ 


h2 


3-V3 
4V2  5 


and 


1-V3 

4V2 


(59) 


with  corresponding  frequency  response  function  given  by 


if©  =  V2  ■  cis(3j©  •  cos2  (j©  ■  cos  (7©  +  i  -J3~  ■  sin  (7© 


(60) 


This  filter  has  a  double  zero  at  £  =  ±Vi. 

Filtering  with  sampling  rate  changes 

Given  an  FIR  filter  E,  the  following  composite  operators  are  defined 
consistently: 


(I2o3*): 


Xt  '  h(:-2n 


(for  analysis ) 


(61a) 
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and 


1 T  2) : 


(for  synthesis ) 


(61b) 


where  the  superscript  *  means  the  adjoint  operator.  Notice  that  each 
composite  operator  in  Equations  6ia-b  is  the  adjoint  of  its  companion. 
Indeed,  one  has 


s°T  2)*  =  (T  2)*  oS*  =  j  2oS* 


(62) 


Also  notice  that,  in  terms  of  the  frequency  response  functions,  one  has 


(|  2  oH*)  x 

a(5)«-4h 

T 

•  X 

X 

+  H 

[5  +  1] 

•  X 

[5+ijl 

\  /  J 

W  2  [ 

.2, 

.2, 

2  j 

2  Jj 

(for  analysis )  (63a) 


and 


(So|  2)  xf  (l)  =  H(l)  -X(2l)  (63b) 


Consider  (for  simplicity)  the  Haar’s  filters  and  fix  the  sequences  a  =  {an}n 
and  c  =  {cn}n.  Then,  the  fundamental  sets  of  sequences  { cp(n) }  and  { (p(fl) } 
are  given  by 


and 


,(2n) 


+  £ 


(2n+l) 


(64) 


1 

12 


(2 n)  _g(2n+l)j 


(65) 


And,  subsequently,  one  can  compute  for  analysis  the  coefficients  produced 
by  the  operations: 
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1.  Averaging: 

(j  2oH*j  a  =  J2^ia2n  +  a2n+i)  e(n)  =  {  (a,  (p(n))  }n  (66) 

and 

2.  Differencing: 

(|2oH*)a  =  ^^-(a2n-a2n+1)e(,,)  =  {(a,cp(?!))  }  (67) 

n  V  2  71 

And,  for  synthesis,  the  orthogonal  series  expansions  produced  by  the 
operations: 

3.  Spreading: 


°T  2)n  —  '  ^n-2£ 

l 


1 

7J'Cm  •  77=2ni 

1 

- •  r* 

,  n=2m+l 


(68) 


(Hot2)c  =  yc„.<p("» 

77 


(69) 


and 

4.  Spreading  and  oscillating: 


n=2m 


n=2m+l 


(Hof  2)  C  =  '  CP 

77 


(n) 


(70) 


(71) 


Applying  operation  1  followed  by  operation  3,  one  obtains  the  process 
(analysis  +  synthesis)  for  the  coarse  details  in  the  digital  signal  a  as 
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(Ho  1 2)o(  1 2o  E*)  a  =  0a,  cpw)  cp(n)  (72) 

77 


and,  applying  operation  2  followed  by  operation  4,  one  obtains  the  process 
(analysis  +  synthesis)  for  the  fine  details  in  the  digital  signal  a  as 

(Ho  1 2)o(  f  2o  H*)  a  =  0a,  cp(n))  cp(n)  (73) 


But,  together,  the  two  sets  in  Equations  72  and  73  form  a  complete  ortho¬ 
normal  family  in  the  space  f2.  Therefore,  the  analytical  expression  for 
splitting  the  digital  signal  a  in  terms  of  coarse-detail  and  fine-detail  com¬ 
ponents,  respectively,  may  be  written  as 

a  =  0a,  <p'">  +  0a,  <(,<">  (74) 

77  77 


Consequently,  the  coarse  details  come  out  of  a  low-pass  filter,  while  the 
fine  details  come  out  of  the  associated  high-pass  filter.  The  down-sampling 
operation  in  Equations  66  and  67  introduces  aliasing.  Subsequently,  the 
up-sampling  operation  in  Equations  68-69  and  70-71  eliminates  that 
aliasing  and,  hence,  allows  the  reconstruction  of  the  original  signal.  This 
decomposition  constitutes  a  fundamental  filter  bank  (Figure  1).  Other 
filter  banks  can  be  constructed  on  the  basis  of  different  filters  and  com¬ 
pound  filters  based  on  this  fundamental  arrangement. 


Figure  1.  Fundamental  filter  bank. 
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Vetterli’s  conditions 

Perfect  reconstruction  of  the  original  digital  signal  can  be  achieved  upon 
successive  analysis  followed  by  synthesis  algorithms  if  the  underlying  filter 
satisfies  the  Vetterli’s  conditions  stipulated  as  follows: 

Let  the  vector  V(^)  be  defined  as  having  two  components.  The  first  com¬ 
ponent  is  the  frequency  response  function  of  the  low-pass  filter  E,  H(Q, 
and  the  second  component  is  the  frequency  response  function  of  the 
associated  high-pass  filter.  The  Vetterli’s  conditions  for  the  filter  are 

||V©f=2  (75a) 

in  other  words,  the  length  of  the  vector  V(^)  =  V2  ,  and 

{V(l),V$+%))  =  0  (75b) 

in  other  words,  V©  orthogonal  to  V(^+V2). 

Scaling  and  wavelet  functions 

Given  a  filter  defined  by  its  coefficient  sequence  h  =  {hm}m,  the  funda¬ 
mental  definition  of  scaling  function  relates  the  function  to  itself  at  two 
contiguous  levels  of  resolution  by  the  expression 


cp(x)  =  V2  • 


XX  cp (2x-m) 


(76) 


in  which  the  filter  satisfies  the  condition 


£ft,„  =  V2 


(77) 


It  can  be  shown  that,  through  the  satisfaction  of  the  Vetterli’s  conditions, 
the  set 


{cp (  x-n); 


oo  <n< 


+ooj 


(78) 


is  an  orthonormal  family  in  the  space  L2(R). 
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Similarly,  the  fundamental  definition  of  wavelet  function  relates  the  func¬ 
tion  to  the  associated  scaling  function  at  two  contiguous  levels  of  resolu¬ 
tion  by  the  expression 


\p(x)  =  V2~  •  -<p(2x-m)  (79) 

m 


The  mother  wavelet  rp(x)  is  compactly  supported  if  the  corresponding 
filter  scaling  function  cp(x)  is  compactly  supported.  In  this  case,  it  can  also 
be  shown  that 


J  rp(x)-  rp(x-s)dx  =  (j  2)  h*(fy  =  e(0)  =  8(s) 


(80) 


In  other  words,  the  set 


|rp(x  — n);  — oo<n<+oo} 


is  an  orthonormal  family  in  L2{R). 
Again,  by  the  same  argument,  one  has 


f  cp(x)  •  rp(x  — s)  dx  1  =  (j,  2 

J"°°  JxeZ 


h*h  =0 


(81) 


(82) 


after  Fourier  transforming  into  the  frequency  domain,  so  that  the  family  of 
scaling  functions  at  a  given  level  of  resolution  (different  translations)  and 
the  family  of  associated  wavelets  at  the  same  level  of  resolution  (different 
translations)  are  mutually  orthogonal  to  each  other.  Moreover,  it  can  be 
shown  that  the  system  of  child  wavelets  {  2k/2  ■  ip(  2kx  -  n):  -oo<  k,n  <+oo} 
[where  k  is  the  level  of  resolution,  n  is  the  translation,  and  the  built-in 
amplitude  2fc/2  is  a  normalization  factor  for  unit  energy]  is  a  complete 
orthonormal  family  in  L2(R ).  Therefore,  this  wavelet  set  constitutes  a  basis 
that  spans  the  space  of  the  finite-energy  functions  in  L2{R). 
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4  Multi-Resolution  Analysis 

Let  Vj  be  the  closed  subspace  of  L2{R )  spanned  by  the  scaling  functions  cp jk, 
for  a  fixed  integer  /,  where  cpy/v)  =  2// 2  •  cp(  2 /  x  -  /c).  The  sequence  of  sub¬ 
spaces  {v; }  defines  a  multi-resolution  analysis  with  the  following 

l  J  J  j=— 00 

properties: 

►  The  subspaces  are  nested: 

■■■cV^1cV0cV1c---  (83) 

►  The  relations  between  coarse  and  fine  scaling  functions  are  determined 
by  the  following  dilation  and  translation  laws,  respectively: 

(p(x)eVj  =>  cp(2x) e Vj+1 ,  jeZ  (84) 

(P(x)gFo  cp(x-n)GVo,  ne  Z  (85) 

►  The  sequences  of  the  Vj  increase  to  all  of  L2(R )  and  decrease  to  the 
empty  space. 

►  The  set  of  integer  translates  {cp(x-A:)}^  ^  is  an  orthonormal  complete 
set,  i.e.,  a  basis  for  y0. 

The  projection  of/[x)  onto  the  subspace  Vj  is  an  approximation  at  the/1' 
level  of  resolution  given  by  f  (x)  =  ^c;1  •  q>jk  ,  with  the  coefficient 

k 

sequence  cjk  given  by  the  correlations 

Cjk  =  ( f’  Tjfc )  =  f  fM  •  T  jk  dx  (86) 

The  mother  wavelet  \p  may  be  derived  from  the  father  scaling  function  cp 
by  the  expression 

iV  ~  k  _ 

\p(x)  =  -j2-^hk  -cp(2x-A:),  where  hk=(- 1)  • h ±_k 
k= 0 


(87) 
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in  which  N  is  sufficiently  large. 

The  orthogonal  complement  of  Vj  in  Vj+i  is  given  by  the  subspace  Wj  in 
such  a  way  that 

Vj+ 1  =  Vj  ©  Wj  ,  j  G  Z  (88) 

where  ©  is  the  direct  sum  between  subspaces.  Similarly  to  Vj,  Wj  is  span¬ 
ned  by  another  orthonormal  complete  family  (i.e.,  basis) 

Wjk(x)  =  2%-W(2  jx-k)  (89) 


where  \p(x)  is  the  mother  wavelet.  This  function  is  orthogonal  to  the 
fundamental  scaling  function  with  respect  to  translation: 


(p(x)-ip(x- 


s) 


dx  =  0 


(90) 


Based  on  the  definition  of  scaling  function,  it  follows  that 


V-  =  V 

j  i 


JEewi+t 

k= 0 


>  J>1 


(91) 


Note  that  for  a  fixed  /,  the  set  {cp jk  }  spans  the  whole  function  space  Vj 
and,  therefore,  any  function  can  be  approximated  as  the  projection 


pAf)=  <92> 

k=— oo 

where  Pj  is  the  projection  operator  onto  the  subspace  Vj . 

Figure  2  shows  a  conceptual  view  of  the  multi-resolution  analysis  pro¬ 
cedure  on  a  given  digital  signal. 
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Figure  2.  Conceptual  multi-resolution  analysis. 
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5  Homogenization  via  Projection  Method 

Let  a  discretized  boundary  value  problem  expressed  in  the  function  space 

Vj+i  be 


T  TT  —  F 

^j+ 1  u  j+ 1  —  rj+ 1 


(93) 


in  which  the  operator  Lj+1  acts  on  the  space  Vj+i  spanned  by  the  basis  func¬ 
tions  at  the  (/+i)st  level  of  resolution.  Consider  the  direct  decomposition  of 
the  Vj+i  space  into  the  scaling -function  and  wavelet-function  spaces: 

vj+1  =  Vj  ©  Wj  (94) 


The  fine  and  coarse  components  of  Vj+i  can  be  separated  using  discrete 
projecting  operators  as 


wj '  Uj+i  = 


Uh 

u  j+ 1 

ue. 

u  J+l 


(95) 


where  Uj+i  is  split  into  a  coarse-detail  component  U'j+1  as  a  projection  onto 

Vj ,  and  a  fine-detail  component  a  projection  onto  Wj ,  and  Wj  is  the 

wavelet  transformation  operator.  In  the  case  of  the  simple  Haar  basis,  for 
example,  the  discrete  form  of  the  operator  Wj  is  given  by 


IV  j  = 


(96) 


where 


1 

1 

0 

0  . 

.  0 

0 

1 

0 

0 

1 

1  . 

.  0 

0 

V2 

0 

0 

0 

0  . 

.  1 

1 

(97a) 


{j2i^1x2 
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and 


1 

-1 

0 

0  . 

.  0 

0 

1 

0 

0 

1 

-1  . 

.  0 

0 

V2 

0 

0 

0 

0  . 

.  1 

-1 

(97b) 


Recognizing  that  operator  wj  is  orthogonal,  the  transformation  of 
Equation  93  renders 


(98) 


(99) 

Applying  static  condensation,  the  reduced  form  of  Equation  98  in  terms  of 
the  coarse-detail  component  gives 

Lj+1  ■  Uj+1  =  Fj+1  —  L2 1  •  L±1  ■  Fj+1  (100a) 

where 

Lj+1  =  ^  L22  —  L21  •  L±1  ■  L12).+i  (100b) 

Then,  the  solution  of  Equation  100a  may  be  decomposed  into  coarse-  and 
fine-detail  components  before  proceeding  recursively  to  the  next  scale,  i.e., 


Wj  ■  Lj .  1  •  Wj  )  I  Wj  ■  Uj , , 


^11 


^21 


^12 


^22 


Uh 

u  j+i 

u(. 

UJ+ 1 


=  WJ  •  Pj+ 1  = 


F 

rj+ 1 

F- 
rj+ 1 


where 


F 

rj+i 

F- 

rj+ 1 


Qj  '  pj+ 1 

Pj  ■  pj+i 
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The  multi-scale  homogenization  procedure  consists  of  generating  and 
solving  this  sequence  of  homogenized  equations  followed  by  the  corre¬ 
sponding  decompositions. 
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6  Example  of  Signal  Decomposition  and 
Synthesis 

Next,  a  numerical  example  of  sequential  analysis  and  synthesis  of  a  simple 
signal  is  presented  (MathSoft  2000).  Figure  3  shows  a  bounded  oscillatory 
signal  with  high-frequency  components  near  the  origin  of  time.  Various 
multi-level  fine  details  in  the  associated  signal  analysis  are  shown  in 
Figure  4.  The  coarse  components  appear  near  the  bottom  of  the  figure,  and 
the  fine  details  are  shown  at  increasing  level  of  resolution  as  the  reader 
moves  up  within  the  figure.  The  reconstructed  original  signal  appears  near 
the  upper  horizontal  margin  of  the  figure.  The  corresponding  sequential 
reconstruction  process  is  shown  in  Figure  5.  Notice  how  the  progressive 
addition  of  details  eventually  recovers  the  whole  signal.  The  resultant 
signal  at  some  intermediate  level  of  resolution  may  be  considered  by  the 
analyst  as  sufficient  for  practical  applications,  thereby  economizing 
representation  and  calculation  efforts. 

The  sequential  exercise  of  decomposition  and  reconstruction  is  repeated  in 
Figures  6-8  for  the  same  signal,  but  with  superimposed  Gaussian  noise. 
Note  that  the  highest-resolution  details  are  essentially  pure  random  noise. 
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EXAMPLE  1.  MULTI-RESOLUTION  ANALYSIS: 

Trun  :=  1.024  AT:=0.001  thisWave  :  =  coiflet(  12) 

N;_Trun  i  :=0..  N—  1  x  := A T  i  N=  1.024103 

AT  1 

doppler(x)  1  -  x)  -sin  2.1- — - - 

x-|-  0.05 


x 

y.  :=doppler  — — 
1  Trun 


MaxD WTLevel( y )  =  1 0  J  :=6  j  :=0„  J-  1 


Figure  3.  Signal  with  initial  high-frequency  oscillations. 


[4]  Reference:E:\Program  Files\MathSoft\Mathcad  8  Professional\HANDBOOK\wavelets\Mra l.mcd 


mra(v,J, filter)  := 


w<—  dwt(v,  J,  filter) 

Z  ,  ,  ,«-0 

M<0>  <— idwt(put_smooth  (Z,  J, get_smooth  (w,J)),J, filter) 
for  qje  J..  1 

z  /  x  ,«-0 

rows(  v)  —  1 

M< J 1  q’>  <—idwt(put_detail(Z,qj,get_detail(w,qj)),J,  filter) 


M 


M  :=mra(y,  J, thisWave  ) 
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S6.  :=Mn  . 

i  0,i 

D6.  :  =  M,  . 

i  i,i 

D3.  !=M.  . 

i  4,i 


D5.  :  =  M  . 

i  2,i 

D2.  :=M,  . 


D4.  :  =  M,  . 
i  3,i 

Dl.  :=M,  . 

i  6,i 


Figure  4.  Multi-level  fine  details  in  signal  analysis.  (Resolution  decreases  from  the  top  to 

the  bottom  of  the  figure.) 
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pf:=l 

mrapprox(y ,  J,  filter) 


Q«—  mra(y,J,  filter) 
for  ie  0..  rows(y)  -  1 

Mo,rQ0,i 

for  je  1..  J 
M.  <-M.  ,  +  Q.  . 

j,i  j—  i,i 

M 


MR  :=mrapprox(y,  J,  this  Wave  ) 


S6.  :=MRn  . 

i  0,1 

S5.  :=MR,  . 

i  1,1 

S4.  :  =  MR.  . 

i  2,1 

S3.  :=MR,  . 

i  3,1 

S2.  :=MR„  . 

i  4,i 

SI.  :=MRt  . 

i  5,i 

SO.  :  =  MR4  . 

i  6,1 

7 

50. -t-6-pf 

-  6 

6-pf 

51.  +  5-pf 

5 

5-pf 

52. -|-4-pf 

4 

4-pf 

S3  .-I-3  pf 
-  •  3 

54. +2-pf 

-  2 

2-pf 

55.  +  lpf 

1 

1  Pf 

56. 

i 


0  200  400  600  800  1000  1200 


Figure  5.  Reconstruction  of  original  signal  in  multi-level  synthesis  process  by 
progressive  increase  of  resolution  through  the  addition  of  fine  details,  starting  by 
the  coarsest  representation  at  the  bottom  of  the  figure. 
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o  :=  0. 1  noisy_doppler  :  =  y-pmorm(rows(y),0,o) 


Figure  6.  Noisy  signal  with  initial  high-frequency  oscillations. 


MN  :=mra(noisy_doppler  ,J,thisWave  )  nd  := noisy  doppler 


S6.  :  =  MN  . 

i  0,1 


D6. 

1 

:=MN,  . 

1,1 
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1 

:=MN  . 

2,1 
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3 ,  i 

D3. 

1 
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4,i 
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1 
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5,1 

Dl. 

1 
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6,1 

nd  -|-7 
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4 
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Figure  7.  Multi-level  fine  details  in  signal  analysis.  (Resolution  decreases  from 
the  top  to  the  bottom  of  the  figure.) 
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MRN  :  =  mrapprox(  noisy  doppler  ,  J,thisWave  ) 


S6. 

1 

:=MRNa  . 

0,1 

S5. 

1 

:=MRN,  . 

1,1 

S4. 

1 

:=MRN0  . 

2,i 

S3.  :=MRN. 
1  : 

S2. 

1 

:=MRN.  . 

4,i 

SI. 

1 

:=MRN,  . 

5,1 

SO. 

1 

:=MRN.  . 

6,i 

50.  -|-6pf 
6-pf 

51 .  -|-5  pf 
5  -pf 

52  -|-4-pf 
4-pf 

53  -|-3-pf 
3  pf 

54  -|-2-pf 

2pf 

55.  -|-  1  pf 

1  pf 

56. 
i 

0 


Figure  8.  Reconstruction  of  original  signal  in  multi-level  synthesis  process  by 
progressive  increase  of  resolution  through  the  addition  of  fine  details,  starting 
by  the  coarsest  representation  at  the  bottom  of  the  figure. 
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7  Simple  Cantilever  with  Heterogeneous 
Harmonic  Elasticity 

In  a  practical  application,  a  simple  cantilever  is  subjected  to  a  concen¬ 
trated  load  at  the  tip.  The  cross  section  of  the  cantilever  is  uniform  along 
its  length,  but  the  material  elasticity  is  highly  heterogeneous.  In  fact, 
Young’s  modulus  of  elasticity  is  modeled  as  a  constant  mean  value  plus  a 
harmonic  fluctuation  about  the  mean  whose  frequency  is  a  controlled 
parameter. 

The  exact  static  displacement  response  of  the  cantilever  was  obtained 
through  the  application  of  the  virtual  work  principle  and  appears  repre¬ 
sented  in  Figure  9.  Figure  10  shows  the  multi-level  fine-detail  decompo¬ 
sition  with  the  coarsest  component  represented  toward  the  bottom  of  the 
figure.  The  reconstructed  original  signal  appears  near  the  upper  horizontal 
margin  of  the  figure.  The  corresponding  sequential  reconstruction  process 
is  shown  in  Figure  11.  Notice  how  the  progressive  addition  of  details  even¬ 
tually  recovers  the  whole  signal.  However,  in  this  case,  interestingly,  all 
levels  of  detail  resolution  have  a  significant  contribution  to  the  signal 
synthesis,  particularly  at  the  boundary  ends,  where  both  the  essential 
boundary  conditions  (displacement  of  the  cantilever  tip  at  the  right  end) 
and  the  natural  boundary  conditions  (reactive  shear  force  and  bending 
moment  enforcing  the  zero  displacement  and  slope  at  the  support)  are 
fundamental  for  the  reconstruction  at  all  stages  of  resolution.  This  is  also 
true  even  when  only  the  coarsest  component  of  displacement  response  is 
considered. 
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EXAMPLE  2.  HETEROGENEOUS  CANTILEVER 

»>  Define  INPUT  variables: 

alpha  0.95  Elastic  modulus  ratio  Eo/Ea 

***  CAREFUL:  alpha=1  INSTABILITY  *** 

n  :=20  Elastic  modulus  wave-number  parameter 

dlambda  Increment  of  response  location  parameter 

128 

this  Wave  :  =  coiflet(  12) 

»>  Compute  RESPONSE  along  cantilever: 

nsubs  :  =  127 

i  :=0..  nsubs 

A  :=i  dlambda  Response  location  parameter 


muO.  := 
1 


1 -I- alpha-cos(2-Jt  -n-i;) 


■d^ 


M-matrix 


mul. 

1 


n. 

i 


Jo 


_ \ _ 

1  +  alpha -cos (2-ji  -n-i;) 


mu2.  := 
1 


1  -i-  alpha -cos  (2-jc  -n-1;) 


-d^ 


5.  :=A.  muO.  -  1-|-A.  mul.-|-mu2. 

Ill  11  1 
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MaxD WTLevel( 8 )  =  1 0  J:  =  6  j:=0..  J-l 


Q  Reference: E:\Program  Files\MathSoft\Mathcad  8  Professional\HANDBOOK\wavelets\Mral.mcd 


mra(v,J,  filter)  := 


w<— dwt(v,  J,  filter) 

Z  ,  ,  ,«-o 

rows(  v)  —  1 

M<0>  «— idwt(  put_smooth  (Z,  J,  get_smooth  (w,J)),J, filter) 
for  qje  J..  1 

z  ,  ,  ,«-0 

rows(  v)  —  1 

M< J  1  tu>  <—  idwt( put_detail  ( Z ,  qj ,  get_detail  ( w,  qj ) ) ,  J,  filter) 


M 


M  :=mra( 8 ,  J,  thisWave  ) 
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<0> 


S6  :=  M 


D6  := 

D3  := 


<1> 

T  <2> 

/  t\<3> 

D5  :=  M 

D4  :=  M 

<4> 

/  T  <5> 

/  T  <6> 

D2  :=  M 

D1  :=  M 

Figure  10.  Multi-level  fine  details  in  displacement-response  analysis.  (Resolution  decreases 
from  the  top  to  the  bottom  of  the  figure.) 
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pf:=l 

mrapprox( y ,  J,  filter)  :=  Q«—  mra(  y ,  J,  filter) 

for  ie  0..  rows(y)  -  1 

Mo,rQ0,i 

for  j  g  1 ..  J 
M.  M.  ,  .  +-  Q.  . 

J  > 1  J  -  1 ,  i  JO 

M 


MR  :=mrapprox(  5 ,  J,  this  Wave  ) 

T  <o>  T  < l >  T  <2>  T  <3> 

S6  :=  MR  S5  :=  MR  S4  :=  MR  S3  :=  MR 

T  <4>  <5>  <6> 

S2  :=  MR  SI  :=  MR  SO  :=  MR 
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i  :=0..  nsubs 


Figure  11.  Reconstruction  of  exact  displacement  response  in  multi-level  synthesis 
process  by  progressive  increase  of  resolution  through  the  addition  of  fine  details, 
starting  by  the  coarsest  representation  at  the  bottom  of  the  figure. 
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Appendix  A:  Summary  of  Reproducing  Kernel 
Representation 


A  uni-dimensional  field  u(x)  may  be  expressed  as  (Liu  and  Chen  1995) 


u(x)  =  uw(x)  +  Jv<lx>[M]_1|^Ju(y)cp 


x-y 


ar 


dy 


x-y 


an 


(Al) 


dy 


where  uw  is  the  part  of  the  field  that  may  be  obtained  by  multi-resolution 
wavelet  reconstruction,  cp(-)  is  the  father  scaling  function,  ao  is  the  dilation 
parameter,  and  the  window  moment  matrix  [ M]  is  given  by 


M  = 


r  1 

I  9 

x-y 

dy 

f  y  -(p 

f*-y| 

mo  mi\_ 

%  j 

ao  , 

m1  m1± 


/y(p 

J  v 


x-y 


an 


dy  f 

J  v 


y2-<\ p 


x-y 


On 


dy 

dy 


(A2) 


The  factor  <1  x>[M]  1|’*’J  inside  the  integrals  in  Equation  Al  may  be 
expanded  in  terms  of  a  correction  function  C(x,  y,  ao),  given  by 


C(x,y,a0)  =  C1  (x)  +  C2  (x)  • 


x-y 


an 


(A3a) 


where 


C1(x)  =  -  —  - 
an 


m 


11 


(A3b) 


2V  ’  a* 


m. 


(A3c) 


and  A  =  det  (M). 
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Expanding  the  function  uw(x)  in  the  wavelet  subspaces,  one  obtains  the 
truncated-series  approximation: 


OO  ^2 

uW  W  =  X)  S  (U>  Wmn )  Wrnn  (*  )  (A4) 

m= 0  n=n1 

where  wavelet  child  \pmn(x)  corresponds  to  resolution  m  and  translate  n. 
Substituting  Equations  A3a  and  A4  into  Equation  Ai,  one  gets  the 
expression 


u(x)  =  fvc(x,yao)-<P  '  u{y)  dy 


00  n2 


+  E  E 

m= 0  n=n± 


f  u{y)wmn(y)dy 


(A5) 


Wmn(X)-Wmn(X) 


where  the  functions  (x)  are  the  approximations  of  the  wavelet  func¬ 
tions  \pmn  (x)  obtained  from  the  reproducing  kernel  reconstruction.  In 

fact,  the  first  term  in  Equation  A5  defines  the  reproducing  kernel  approxi¬ 
mation  (low-frequency  components),  whereas  the  second  term  in  that 
equation  gives  the  contribution  of  the  high-frequency  components.  [For 
the  low-frequency  components,  ipmn  (x)  ~  (x) .] 
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